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Abstract 

A mutator is an allele that increases the mutation rate throughout the genome by dis- 
rupting some aspect of DNA replication or repair. Mutators that increase the mutation 
rate by the order of 100 fold have been observed to spontaneously emerge and achieve 
high frequencies in natural populations and in long-term laboratory evolution experi- 
ments with E. coli. In principle, the fixation of mutator alleles is limited by (i) competi- 
tion with mutations in wild-type backgrounds, (ii) additional deleterious mutational load, 
and (iii) random genetic drift. Using a multiple locus model and employing both simula- 
tion and analytic methods, we investigate the effects of these three factors on the fixation 
probability Pfi^ of an initially rare mutator as a function of population size A^, benefi- 
cial and deleterious mutation rates, and the strength of mutations s. Our diffusion based 
approximation for Pfi^ successfully captures effects (ii) and (iii) when selection is fast 
compared to mutation {jj,/ s <^ 1). This enables us to predict the conditions under which 
mutators will be evolutionarily favored. Surprisingly, our simulations show that effect (i) 
is typically small for strong-effect mutators. Our results agree semi-quantitatively with 
existing laboratory evolution experiments and suggest future experimental directions. 
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The most evolutionarily important characteristic that an individual inherits from its parents is 
the average number of offspring that it will leave in the next generation, i.e. its fitness. But, is 
fitness the only evolutionarily relevant heritable trait? The ultimate fate of an individual depends 
not only on its immediate properties, but on those of its entire lineage of descendants. Therefore, 
the genetic system that shapes the statistical properties of this lineage is also an evolutionarily 
relevant, selectable trait. 

In this article we study one such property, namely a globally elevated mutation rate. In practice 
this property is inherited via a mutated copy of a gene, called a mutator allele, involved in DNA 
copy or repair. We ask the following basic question: What is the fixation probability of an initially 
rare mutator? This is a generalization of the classic pop ulation genetic calculation for the fixation 



probability of a static mutant with selection coefficient s (|FisherI . 



19301). If the fixation probability 



of a mutator allele differs from that of a neutral one (i.e. 1 /N), then the average mutation rate of 
the population will be under selective pressure. 

The selective forces acting on mutators is not purely a t heoretical issue- 
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tissues of multicellular sexual organisms comprise populations of asexually reproducing cells pos- 
sessing opportunities for an increased gr owth rate. Correspondingly, tumoregenesis has been asso- 
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2OO2I : 



Mao et ali 



196C)have 



GiRAUD et al. ■ 



19971) indicate 



that mutators achieve fixation because of the adaptive mutations they generate and not because of 
any intrinsic fitness advantage. Thus, selection on mutator alleles occurs via an indirect mecha- 
nism. One of the goals of our work is to make semi-quantitative contact between our model of 
indirect selection and the existing data of mutator fixation in laboratory experiments. 

The evolution of mutation rate is a problem that dates back to the 1930's. The general issue was 
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articu 



(19671) and 



ated by 



Leigh Jr 



Sturtevant (|1937|) . and important theoretical contributions date back to 



KiMURA 



1970I). Theoretical studies proli ferated during the last decade, and the field 

Denamur and MATig ((2006) . Given the 



is reviewed by ISniegowski et all iQOOQ) and also by 
abundance of existing theoretical articles, it is critical to understand how our work relates to and 
improves upon this body of literature. We address this issue in detail in the Discussion section. For 
now, we merely provide a brief sketch. First, we neglect the complicating influences of recombina- 
tion and environmental fluctuations. This allows for a direct and comparatively precise treatment 
of the simplest situation: a strictly asexual population adapting in a constant environment. Even 
this simplest scenario has rich and often counterintuitive behavior. Secondly, our methods natu- 
rally treat both strong (e.g. 100 fold) mutators and weak modifiers of mutation rate. Thirdly, unlike 



most previous work, we combine fully stochastic simulations with an ana 
alytic results for weak modifiers are a generalization of previous work by 



ytic approach. Our an 



Andre and Godelle 



(|2006h . but we find that both approaches often fail to match simulations. However, our work for 
strong mutators does match simulations over the expected parameter range. The simulations thus 
provide vital checks and guidance for the analytic approach. Conversely, the analytic approach 
deepens our understanding of mutator fixation and makes predictions in parameter regimes that 
are computationally inaccessible via simulation. Finally, unlike previous work, our diffusion based 
analytic approach captures the effects of random genetic drift. This not only allows for exploration 
of regimes where random drift is important, but also a quantitative understanding of when it can be 
neglected. 

The outline of this article is as follows. We begin with a heuristic discussion of mutator dynam- 
ics. Next, we construct and simulate a stochastic model of asexual populations that include mutator 
alleles. We do not explicitly allow for the formation of mutators, merely the competition between 
mutators and wild-type strains once mutators arise. Afterward, steered by the outcome of simula- 
tions, we develop a quantitative understanding of the results of the stochastic simulations. Although 
a full mathematical treatment turns out to be intractable, we are able to devise an approximation 
scheme that captures many features of the simulation results. We then solve our approximation 



scheme, both numerically and analytically. The resulting expressions allow a comparison to the E. 



5 



Table 1: Commonly used notation. 



Symbol Usage 

N Total population size 

Wild-type mutation rate per genome 

1^1+ Mutator mutation rate per genome 

U Mutation rate into mutator state 

L Length of genome 

b Number of I's in genome 

5 Fraction of mutations that are lethal 

X Mutator frequency 

p. = (1 — x)ij,^ + Average mutation rate per genome 

R= ji^/ jj.^ Mutator strength 

r = b/L Growth rate per individual per simulation time-step 

s = 1/b Selection coefficient of non-lethal mutation 

a = I — b/ L Fraction of O's in the genome 

Qfe = a(l — (5) Fraction of mutations that are beneficial 



coli experiments of Lenski and co-workers (|SnieG0WSKI et ali 



19971). 



HEURISTIC ANALYSIS 

In this section, we briefly explain the conceptual factors underlying mutator fixation. The equa- 
tions in this section should be considered merely as heuristic guides and not formal results. 

Since mutator alleles do not directly affect fitness, their dynamics must be guided by association 
with other genes which do have a direct fitness effect. In asexuals, all loci sha ring the same genome 



with a sweeping beneficial mutation will also become fixed via "hitchhiking" (|Maynard-Smith and HaighI . 



19741) . Whereas most alleles hitchhike completely passively, the mutator allele plays a somewhat 
active role in facilitating its own hitchhiking by increasing the probability of a beneficial mutation 
elsewhere in the genome. This well known mechanism occurs in our simulations and is evident in 
Figlll 

At the same time, the wild-type subpopulation also generates advantageous mutations. When 
this occurs, mutators become extinct due to fixation of their counterpart wild-type alleles. Although 
the wild-type generates mutations more slowly on a per capita basis, if it vastly outnumbers the 
mutator subpopulation, then the total mutation rate in the wild-type background may be larger. 



Along these lines, it is tempting to think of the number of mutators as initially constant, and that 
the mutator will achieve fixation if and only if it generates a sweeping beneficial mutation before 
the wild-type background does. This means that 

p _ _ XojJ,+ 

Jl Xo/i+ + (1 - Xo)/i- 

where Xo is the initial frequency of mutators and yU+ (/i_) is the genome-wide mutator (wild-type) 
mutation rate. This equation has striking qualities. First, it is independent of the following prima 
facie important parameters: population size N, selection coefficient of mutations s, and the fraction 
of mutations which are beneficial versus deleterious. Secondly, and more subtly, the equation is 
explicitly frequency dependent. It will turn out that EqU] arises as a limiting form of our analytic 
expression, but does not typically match the results of simulations. 



In contrast to the frequency dependent EqlH a classic result from population genetics (|FisherI . 



19301) is the fixation probability of a mutant with a simple selective advantage: 



This result holds for haploid populations using Moran process dynamics, and merely requires fac- 
tors of two in the exponents to handle diploids or Wright-Fisher dynamics. In Eql2l Pfix depends 
on the frequency of mutants only via the product Nxo, i.e. the initial number of mutants. Thus, 
Eqs Jll2l scale differently with population size. The form of Eql2l implies that (when NS ^ 1), 
Pfix ~ 1 — e~^^°^ ~ 1 — (1 — S)'^'^° and we can think of each mutant as an independent "trial" 
with fixation probability S. In other words, if the fraction Xo is kept constant and N is increased, 
Eqd] says that Pfi^ should remain unchanged whereas Eql2] says that Pjj^ should increase. On the 
other hand, if Nxo is held constant as is increased, EqH] predicts a decrease in Pfi^ whereas 
Eql2l predicts that Pfix remains unchanged. Since mutators achieve fixation by hitchhiking with 
mutations which are themselves governed by Eql2l perhaps we should a priori view EqU] with sus- 
picion. Indeed, our simulation data and analytic methods will show that mutator fixation is often 
governed by an equation with the form of Eqj2l 
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While Eqdl completely neglects deleterious mutations, they are the basis for another heuristic 
line of thought. In any realistic biological population, regardless of how maladapted, deleterious 
mutations vastly outnumber advantageous ones. Because of this, upon first thought, one might 
think that the mutator allele will do more harm than good and therefore be selected against. Al- 
though it is true that an elevated mutation rate will quite likely cause an immediate decrease in 
the population's mean fitness, evolution does not always act to maximize this quantity. The situ- 
ation is understood more clearly in the following game theoretical context. A beneficial mutation 
often greatly increases the probability that a lineage will achieve complete evolutionary success 
by sweeping through the entire population, whereas a deleterious mutation only slightly decreases 
the low probability of a neutral sweep. More quantitatively, we can think of the "payoff" for a 
sweeping advantageous mutant as the entire population size N. For this to occur, the mutator must 
generate a beneficial mutation which must then survive in spite of random drift. In contrast, the 
payoff for a deleterious mutant is merely a single individual who is destined to die out with near 
certainty. The mutation strategy is favored when its expected payoff is greater than zero, i.e. 

N ■ 7r{s) ■ Hben - 1 ■ l^del > (3) 

where 7r(s) is the fixation probability of a simple mutant, given by Eql2]and Hhen il^dei) are the 
beneficial and deleterious mutation rates, respectively. Note that this expression weights beneficial 
mutations ■ 7r(s) times more heavily than deleterious ones, underscoring their asymmetric effects. 
Later in this article, we show that Eql3]also follows from a rigorous mathematical analysis. 

Thus far we have argued that the fate of mutators is in principle limited both by competition 
with wild-type and by their increased load of deleterious mutants. Additionally, random genetic 
drift is commonly a potent force acting on rare subpopulations. Each mutator begins its existence 
selectively neutral. It can be shown that random drift eliminates neutral alleles from the population 
with a high probability = 1 — and tha t the a verage time taken to do so is merely ~ In(A^) 



generations (|Crow. J.F. AND KiMURA. M, 



19701) . Although we cannot write down a "back of 



the envelope" estimate of this effect, we will later derive a formula that fully incorporates random 
drift and specifies the parameter regimes in which it dominates mutator fixation. 



Our analytic work results in a formula for the mutator fixation probability in terms of simple 
parameters. Examining this expression yields a quantitative sense of the relative importance of 
random drift, deleterious mutations, and beneficial mutations. This allows us to define "strong- 
effect" and "weak-effect" mutator regimes in terms of the model parameters. In the strong-effect 
regime, mutations in the wild-type background do not affect mutator success and our analytic ap- 
proach works well. In the weak-effect regime, mutations in wild-type backgrounds are predicted 
to be the dominant influence on mutator fixation. However, in the case of weak-effect mut ators, 
we will show that our analytic approach, like existing work by I ANDRE and GodelleI (|2006|) . typ- 
ically overestimates the competitive effects of mutations in wild-type backgrounds. When this is 
true, EqUlprovides a poor description of mutator fixation. We now turn toward a discussion of our 
stochastic simulations, that provide an invaluable reference to which we compare our analytic work. 



DESCRIPTION OF STOCHASTIC SIMULATIONS 

We model haploid asexual populations of fixed size N undergoing stochastic processes of birth, 
death, and mutation. Initially, a fraction Xo ^ 1 of the population are mutators and all individuals 
have the same fitness. The birth-death-mutation process is iterated until the population consists 
entirely of either mutators or wild-type. Transitions between the mutator and wild-type states are 
not allowed. We do not model environmental changes explicitly, thereby assuming that the process 
of mutator fixation occurs on a time-scale much shorter than that associated with environmental 
changes. 



19921) . The 



Our stochastic simulations are based on the well known "Moran Process" (|MORANl . 
following sequence of actions occurs every discrete timestep: 

1 . A randomly selected individual is chosen as a potential parent. 

2. The chosen individual gives birth with probability proportional to its fitness. If it does not 
give birth, the simulation advances to the next timestep. 

3. A randomly chosen individual, other than the baby, is killed. 
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4. The baby undergoes a deleterious (beneficial) mutation with probability equal to its deleteri- 
ous (beneficial) mutation rate. This mutation rate of course depends on whether the baby is a 
mutator or a wild-type. Mutations between mutator and wild-type alleles are not allowed. In 
effect, this assumes that mutators are generated on a time-scale much longer than that of the 
entire "competition experiment". 



W e model the genorne of ea ch individual as a string of L bits (|Crosby 



19961 : 



TSIMRING et al. ■ 



19701 : 



Woodcock and Higgs 



19961). A fraction 6 of these bits correspond to critical sites in the genome 



that, when mutated, cause a lethal phenotype. In this case, the baby is never born, and the sim- 
ulation simply advances to the next time-step. Changing the value of 5 in effect allows for some 
adjustment of the distribution of deleterious mutational effects. The birth probability per unit time, 
which we denote r, is proportional to the log-fitness of the chosen individual and equals the fraction 
of I's in the genome, denoted by b/L. Key parameters are a = 1 — b/L and «£ = (!— 6)a, i.e. the 
fraction of sites that would be beneficial if mutated. Thus, all non-lethal mutations have the same 
strength and genes do not interact. This scheme for assigning fitness to genotypes is known as a 
"multiplicatiye Fujiy ama" fitness landscape, and is the K = version of Kaufman's "NK" model 



(|KauffmanI . 



19931). This toy landsca pe is obviously a u s eful mathematical simplifi cation. Ad 



ditionally, recent experimental work by 



Hegreness et al. 



( 2006 ) 



Desai et al. 



(l2007h shows that 



some dynamics of real bacteria and yeast populations can be captured by considering mutations of 
only a single strength. 

Mutation is implemented by "flipping" bits with a probability ^ per bit per birth event, de- 
pending on whether the baby is a mutator (+) or a wild-type (— ). The total number of flips is 
determined by drawing a binomially distributed random number with success probability ^ and 
number of trials L. Each mutation has a probability 6 of being lethal. If no mutations are lethal, the 
number that are beneficial is determined by drawing another binomially distributed random number 
with success probability a and number of trials equal to the number of flips. Unless fi± is 0(1), 
the probability of more than one mutation occurring during a single birth event is negligible and 
we will refer to the genome-wide mutation rate as n±. 

Another useful parameter is s = — a) which, like a and ctg, changes throughout the 
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Figure 1: Some sample runs from simulations where the wild-type mutation rate is zero. The top 
panels depict the number of mutators in the population vs. ^, where r is the birth probability per 
time step which is proportional to the (initial) mean population fitness. The bottom panels show 
the average number of beneficial mutations in the mutator subpopulation. The dark lines resulted 
in fixation of the mutator allele, whereas the lighter lines resulted in its loss. When the mutation 
rate of the mutators is not too large, the mutator hitchhikes to fixation with a single beneficial 
mutation (left panels). When /i+ is larger, many beneficial mutations occur during the fixation 
process (center and right panels). Our analytic approximation scheme assumes that the fixation 
process is triggered by merely the first beneficial mutation to survive drift. Note that in each 
case the population is always far from the fitness maximum when the mutator achieves fixation 
since there are 80 possible beneficial mutations. Parameters are N = 10^, Xq = .005,5 = 0, 
wild-type mutation rate ji^ = 0, and /i+ = 10~^ (left), /i+ = 10^^ (center), /!+ = 1 (right). 
a = A,s = 1/120 (initial values). 
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Figure 2: Averaged results of simulations, and the utility of S^^ as the measure of mutator success. 
When Pfix <^ 1, Pfix increases linearly with Xo (data not shown). The left panels show the (least 
squares) slope of said linear increase when the population is well adapted (bottom) and poorly 
adapted (top) to its environment. The data on the bottom row are quite noisy because of the small 
number of trials resulting in fixation. The panels on the right express the same data, but in terms 
of the effective selection coefficient of the mutator allele obtained by inverting EqlH Whereas 
the values from the left obviously depend on N, the values on the right panels are independent ofN 
when NS^ ^ 1. This suggests that S^, which exposes an underlying simplicity to the simulation 
results, is a more natural measure of mutator success than Pfix- Notice that when the mutator is 
favored, 5^ is always less than the selective advantage s of a single beneficial mutation; this is 
due both to deleterious mutations and loss due to random drift. Parameters are s = 1/120, fX- = 
0, 5 = 0, a = .4 (top) and .008 (bottom). See Supplementary Information for details concerning 
averaging. 
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Figure 3; Dependence on the underlying selective advantage s. The data corresponding to two values 
of s, i.e. two values of L, approximately collapse onto a single curve when and are each 
scaled by s. The scaling of the independent variable underscores the fact that mutator success for 
fixed a is largely controlled by the ratio of timescales for mutation (1 / and selection (1 / s). In 
particular, the sharp decrease in S*^ at large /i+ occurs when these timescales become comparable, 
i.e. when deleterious mutations accumulate in an expanding lineage before it has sufficient time to 
achieve fixation. Parameters are N = 5000, /i_ = 0, a = .4, 5 = 0. 
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simulation as the population evolves. We emphasize that this fitness dependent value of s does 
not represent an epistatic effect. Rather, it is a consequence of mutations which result in a fixed, 
additive increment in "log-fitness." 

A consequence of our genomic model is that both the beneficial and deleterious mutation rates 
will be larger than values encountered in biological populations unless L is extremely large. While 
this may seem like an unnecessary and undesirable restriction, it will turn out that our analytic 
results, which readily handle arbitrary values of the mutation rates, are insensitive to these details 
of our bit string simulation model. 



SIMULATION RESULTS 

To simplify matters, we first investigate the case where the wild-type mutation rate is zero; 
results for the more general case will be given later. Fig[T] shows typical runs for this case. These 
graphs make it clear that if the mutator mutation rate, is sufficiently small, the mutator allele 
hitchhikes to fixation with a single beneficial mutation. This simple observation reminds us that 
mutator fixation or loss is not the result of winning the race up the fitness landscape, but rather 
hitchhiking with beneficial mutations. Thus, mutator allele s are better thought of as co nsequences 



Sniegowski etali 



200nh . When /i+ 



of asexual evolution than causes of more rapid evolution ([; 
is larger, the dynamics are more complex. Despite this complexity, we will later show, via the 
success of our analytic approximation scheme, that the fixation process is triggered mostly by the 
first beneficial mutation to escape random drift. 

Dependence on Fig j2] presents simulation results for three different population sizes and 
two different degrees of adaptation. The fundamental measured quantity is the fixation probability 
Pfix of an initially rare mutator. When Pji^. <^ 1, the mutators are completely independent of one 
another and Pfix increases linearly with Xo (data not shown). To normalize against the effect of Xo, 
we consider the slope of said linear increase, dPfi^/dxa, which equals the mean number of mutator 
descendants left by each mutator, as our preliminary measure of mutator success. Figl2](left panel) 
shows how dPfix/dxo depends on The small and large /i+ limits make qualitative sense: as 
0, the mutator phenotype is "turned off" and therefore neutral, resulting in dPfi^/dXo 1. 

14 



On the other hand when /i+ > 1, a mutation occurs nearly every birth event and the fitness of 
an evolutionary line of individuals takes a biased random walk toward the much lower fitness 
of a completely random genome. Thus, although it is computationally prohibitive to measure a 
negligible fixation probability, it is clear that the mutator allele is nearly lethal at sufficiently large 

Dependence on N, and mutator effective selection coefficient: Figl2]also shows that dPfi^/dx, 
increases with increasing A^. This behavior is incompatible with Eqdl which is independent of N, 
but is fully consistent with Eql2l 

p _ ^ 

We now quantitatively consider whether Eq|2l which applies to mutants with a direct fitness ad- 
vantage, also describes mutators with indirect fitness effects. For this to be the case, the fixation 
probability measured from simulations with differing values of N and Xo would all correspond 
to a single value of S^{a, s, yU.+ , 5). Using the values of Pji^ measured from simulations, we 
used a computer to invert Eql2l thereby obtaining corresponding values of 5^. Figl2l^right) shows 
that, when NS^ ^ 1, there indeed exists an underlying quantity S*^, which we call the "effective 
mutator selection coefficient," that remains invariant as N, Xo, and Pf i^ change. 

There are several advantages to using as the measure of mutator success. First, it allows Eqj2] 
to determine in advance how Pfi^; depends on and Xo, thereby reducing our number of parameters 
by two. Secondly, it allows us to apply aspects of our conceptual understanding of direct mutants to 
the fixation of indirect mutators. For example, when NS^ ^ 1, Pfix for a single mutator becomes 
independent of N, i.e. the notion of a frequency independent per capita fixation probability makes 
sense. Thirdly, the existence of S^, in the sense of Eqj2l invites future questions. For example, 
one may wonder whether S*^, in addition to determining Pfi^, also describes the average dynamical 
behavior of the mutator subpopoulation, e.g. whether ~ e'^'"* when rare. In this article we 

do not apply such an interpretation on S^. Rather, we merely interpret it as a succinct descriptor of 
mutator success. 

Dependence on strength of mutations: Figl3] shows how S*^ depends on the strength of the 
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mutations on our fitness landscape, as measured by s. Figl3](left) shows that as s is increased, 5*^ 
also increases, and reaches its maximum value at a faster mutation rate. FigjStright) demonstrates 
that the curves in the left panel are not as different as they appear: when 5*^^ and yU+ are each scaled 
by s, the curves become nearly identical. This means that is directly proportional to s, and that 
is governed by the single composite parameter rather than /i+ and s separately. Thus, an 
examination of the simulation data has allowed us to reduce our number of parameters by three. 

INSTANTANEOUS SINGLE LOCUS APPROXIMATION (ISLA) 

Stochastic simulations provide valuable signposts along the way to understanding mutator fixa- 
tion. However, a deeper understanding, as well as the ability to probe computationally prohibitive 
regions of parameter space, requires an analytic approach as well. At a given time, the state of 
the population is fully specified by (i) the number of mutators, (ii) the fitness distribution of the 
wild-type subpopulation, and (iii) the fitness distribution of the mutator subpopulation. A complete 
solution to the stochastic process requires an enumeration of the transition probabilities between 
each of these states at each point in time. The problem with such an approach is the extremely 
large number of possible fitness distributions and the correspondingly high dimensionality of the 
resulting governing differential equations. In order to make progress, we note the heuristic rule 
that deleterious mutations are rapidly removed from the population, whereas beneficial mutations, 
and all loci linked to them, become rapidly fixed. This observation motivates the following ap- 
proximations that handle mutations, which are the ultimate source of the aforementioned daunting 
multiplicity of fitness distributions. 

Approximation 1 (Al ): We assume that when a beneficial mutation arises, it instantly becomes fixed 
with a probability given by the classical fixation probability tt of a beneficial mutation in a static, 
homogeneous environment. For our Moran process dynamics, this probability is simply s if s <C 1 
and Ns ^ 1. All loci in the genome in which the beneficial mutation arose also achieve fixation 
via hitchhiking. This represents the most common process by which the mutator allele achieves 
fixation or loss. 

Approximation 2 (A2): The remaining fraction 1 — s of beneficial mutations are simply ignored 
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and treated as if no mutation occurred. This approximation is necessarily somewhat awkward. On 
the one hand, A2 is unnatural in that it allows lineages which are destined to be extinguished by 
random drift to remain in the population and potentially generate their own beneficial mutants. An 
alternative, which we call A2*, is to immediately kill the beneficial mutants which do not sweep, 
which is clearly too harsh. These two alternatives lead to a trivial difference in our formulas, and 
are discussed in Supplementary Material. 

Approximation 3 (A3): Deleterious mutations are treated as effectively lethal, since their descen- 
dants are quickly removed from the population. This results in an effective reduction in the birthrate 
of the mutator strain. 

Since these approximations preclude fitness polymorphism over finite time intervals, they allow 
us to describe the dynamics of the entire population with the single time dependent random variable 
X, i.e. the frequency of the mutator locus. Approximating a; as a continuous variable, and expressing 
time in "generations," the diffusion equation governing P(x, t) is (see Supplementary Information 
for a detailed derivation) 



Each of the three lines in Eq|4]has a straightforward physical interpretation. The first line represents 
"random genetic drift." The second line represents the mutational load of the mutator. The final line 
represents the "decay" of probability from the open interval x E (0, 1) due to beneficial mutations 
that instantaneously sweep. 

An approximation to a limited version of Eqj4]is solved in Supplementary Information. How- 
ever, we can write an equivalent "backward Kolmogorov" equation which is often more mathemat- 
ically convenient than Eq. |4l Defining G{xo, t) as the probability that the mutator has been lost by 



dP 
'dt 




— Na^s + (1 — P 



(4) 
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time t, given that a; = at t = 0, we find 

dG{xo, t) 1 

- -Xo{l - Xo)—G{Xo,t) 



dt N "'dx 



o 



d 

- [1 - ae(l - s)] Xo{l - Xo)-^G{Xo, t) 
Nfj,^aeSXoG{xo, t) + Nfj,_aes{l — Xo)(l — G{xo, t)) (5) 



The backward equation is primarily useful in its steady state form. Defining G{xo,t oo) = 
Goo{xo) and taking the continuum limit, we obtain the ODE 

- l^G 

" Ndxl 

- - [1 - ae(l - S)] ^Goo 

— Nfi+ae-s — h Nn_aeS- — (6) 

1 Xq Xq 

Solution and Analysis Without Wild- Type Mutations: We return for now to the simpler case 
fi_ = 0, deferring until lat er the more general situation. Eq|6] can be solved exactly in terms 



of the Whittaker M function dABRAMOWlTZ and StegunI . 



19651). This exact solution is however 



not immediately instructive (and in any case cannot be generalized to the case of finite wild-type 
mutation rate). It is simpler in practice to solve Eql6]numerically (see Supplementary Information). 
It is also possible to extract some useful information directly from the differential equation. 

First, we note that a simple analysis reveals when the mutator allele will be favored. For nota- 
tional convenience we define the constants 

B = fi+ [1 - ae{l - s)] 
G = fi+a^s 

According to ISLA, the mutator is neutral /or all /i+ when Gooixo) = 1 — Xq. Plugging this into 
Eql6l we find that this requires B = NG, or 

1 1 

" ~ 1 + {N-I)s^ 1 + Ns ^ ^ 
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Figure 4: Behavior near the transition from favored to disfavored mutators. When a^. is greater than 
a critical value ag*"**, the mutator allele is favored (S*^ > 0) for small enough Our analytic 
approach (ISLA) predicts that the transition occurs at [N s + 1 = 1, which agrees extremely 
well with simulation data. Parameters are = 5000, s = 1/120, yU_ = 0, 5 = 0. The number of 
available beneficial mutations are, in order of decreasing mutator success: 10, 5, 3, and 1. 
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Figure 5: Comparison of simulation, numerical solution of Eql6l and the analytic approximation 
EqlTOl The exact numeric solutions to our ISLA Eql6] for different converge to the analytic 
approximation EqlTOlwhen NSfj_ ^ 1 (left). Solutions to Eql6]show, in agreement with simulation, 
that Sf^/s depends on rather than and s seperately (right). Parameters are those used in 
Figs, m [31 

It is easy to check that this expression also holds for the fi^ > case. First note that Eq|7] is 
equivalent to our heuristic guess, EqjBl if Ns ^ 1. Examining EqJTl we see that conditions which 
favor the emergence of mutators (at least when the resident mutation rate /i_ is negligibly small) 
are large population size, potent mutations, and a relatively large fraction of sites that would 
be beneficial if mutated, perhaps due to an environment to which the organism is not well adapted. 
The fact that large tte favors mutators is obvious. The dependence on N is simply a result of the fact 
that as population size increases, the neutral fixation probability becomes an easier benchmark 
to exceed. The qualitative dependence on s is also straightforward in hindsight, given A1-A3: 
increasing s increases the fraction of beneficial mutations that achieve fixation, but does not affect 
the fate of deleterious mutations, all of which are treated as lethal. Also notice that for sufficiently 
large the mutator is always favored, although its fixation probability may be very small: it is 
favored only in the sense that it fares better than a neutral allele whose fixation probability is 1/A^. 
Fig m demonstrates the success of Eq|7]when /i+/s ^ 1. The failure of ISLA for larger /i+/s will 
be discussed later. We next develop approximate solutions to Eql6l with = 0. 
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Strongly Favored Mutators ( NSf^ 1): In this regime, we expect Pfi^ to increase rapidly with 
Xo- Therefore, we expect the loss probability Goo{xo) to decrease rapidly, and 1/(1 — Xo) to differ 
significantly from 1 only when Goo ~ 0. Then, for Xo ^ 1, we can approximately take 1 — Xo ^ 1, 
and the solution to Eql6]with /i_ = is simply 

Gooixo) = e-^^^° (8) 
VB^ + AC-B 
' " 2 

Our approximation is self consistent if indeed Goo decays rapidly, i.e. Nz ^ 1. This solution does 
not satisfy the boundary condition at Xo = 1 since our solution is only valid for Xq ^ 1. Beyond 
this region the structure of the solution is more complicated, which need not concern us here since 
fixation is essentially total in this regime. We then have for the fixation probability of the mutator 

P^x(xo) = l-e-^^^° (iV^»l) (9) 

A comparison with Eql2] shows that, according to A1-A3 and in the limit Nz S> 1, the mutator 
effectively behaves like a simple advantageous mutant with a well defined selection coefficient 

= z: 



_ VB^ + AC-B /i+ 



a/(1 - UeY + 4aeS//iH 



Up, 



iV^„>l (10) 



A comparison of the stochastic simulation data with both a numerical solution of Eql6] and this 
approximate analytic expression (EqfTOl) is given in FigjSl We see that our approximate S^/ s only 
depends on rather than /i and s separately, as we noted in the Simulation Results section. 

For small ^ agS, C ^ B^ and ~ ^/C = ^^^a^s, and thus only advantageous 
mutations are relevant to mutator success. This result (which is directly supported by Figj7]to be 
discussed later), shows that in this regime, random drift, and not deleterious mutations, is the only 
check on mutator success. 

Li the complementary regime where /i+ ^ OgS, \S^\ approaches its maximum value S"* with 
respect to Here, the solution is the same as if the second derivative term, which represents 

21 



random drift, were dropped from Eql6] (see below). Therefore, random drift is irrelevant in this 
regime and deleterious mutations alone limit mutator success, giving 

^ B 1 — ae 

The factor in EqdT] multiplying s is the ratio of beneficial mutations to deleterious and lethal 
mutations. In real biological populations, this ratio is certainly less than one, and hence S** ^ s. 

Marginal Mutators (NS^,, < Ij: We can readily make progress in this regime if Nfi^ ^ 1 and 
N'^jj.j^aeS ^ 1. In this case, the B and C terms dominate Eqj6]and the solution for G^o is simply 

G^{xo)^{l-Xof'*- (12) 

with a fixation probability Pfix{xo) ~ NxqS*. In obtaining this solution, we dropped the second 
derivative term in Eqj6l which could in principle introduce large errors near Xo = 1, where G"{xo) 
from Eqdllis in fact large. Nonetheless, it turns out that Eq [12] satisfies the boundary condition at 
Xo = 1 and thus remains a valid leading order approximation for all Xo- Since Pfi^ is comparable 
to in the present marginal case, we cannot interpret S** as a mutator selection coefficient here. 
Rather, we have Pji^ = + NSf,/2), from which we obtain iVS*^ = 2 °"^^^+^^"^ independent 
of The numerator of this expression makes clear the agreement with our previous estimate for 
the critical value of ae given by EqjTl 

The case where A^/i+ < 1 and NS^ < 1 requires a more lengthy analysis, and is presented in 
Supplementary Information. 

EFFECT OF WILD-TYPE MUTATIONS 

We now turn our attention to the more complicated case when mutations in wild-type back- 
grounds are allowed, i.e. > 0. We begin by solving Eql6]for /i_ > in the large Njj,± limit, 
where the second derivative term can be neglected. Working in this limit simplifies the mathe- 
matics, and is sufficient for illustrating the points that we intend to make. An approximation that 
incorporates the second derivative term and random drift is included in Supplementary Information. 



22 



G-ON=10".n_=0 
&^N=10'. n^/(l_=100 
N=5000. |l_=0 

N=5ooo, =100 

x-^ N=5000. |l^/(l_=IO 
N=5000. |l m =2 




0.0001 



o simulation, N = 1000 
— ISLA, N = 1000 
□ simulation, N = 5000 
ISLA, N = 5000 

^ simulation, N = 10^ 

■ ISLA, N= lo' 




0.01 



100 



Figure 6; Small effect of mutations arising in wild-type backgrounds. ISLA predicts that these mu- 
tations will become important in the weak-effect mutator regime defined by ^ 1^ where 
R = However, the simulation data show that mutations in wild-type backgrounds some- 
times have a negligible impact even in the weak-effect mutator regime. In the panel on the right, 
^Na°'s^ values 18, 3.6, and .18, respectively, as N is increased. Accordingly, ISLA predicts 
a decrease in S*^, but did not change in simulations. The panel on the left shows that benefi- 
cial mutations in wild-type backgrounds eventually decrease for large enough R, though the 
decrease here is smaller than what ISLA predicts. Parameters are a = .4, s = 1/120, 5 = 0, and 
12+/ = 100 (right). 
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In the large limit, 







AT 



I — X, 



oo 



+ Nfji^aeS 




This first order, linear ODE can be solved by standard methods. Defining R = ji^/ ji^, we obtain 



The prefactor in Eq[T3]is identical to our previous expression for the = case (Eqs ll 1112)) 
when Xo <^ 1. Recall that the sign of the quantity ae{Ns + 1) — !^ NaeS — 1 determines whether 
mutators are favored (EqjT]). Therefore, mutations in wild-type backgrounds decrease Pfix when 
mutators are favored and increase Pfix when they are disfavored. This latter effect occurs because 
mutating is generally a losing strategy when ae{Ns + 1) — 1 < (see EqjS]): the small persistent 
cost of deleterious mutations exceeds the huge occasional benefit of a selective sweep. Thus, in 
this regime the wild-type aids the mutator by participating in this losing strategy. 

EqfT3] also determines when R is sufficiently large to ignore mutations in wild-type back- 
grounds. In other words, EqUS] allows us to define natural "strong-effect" and "weak-effect" 
mutator regimes. For weak-effect mutators, ""^^'^^^"^j""'' ~ Na^s ^ R, and Eq{T3] reduces to 
Pfix = XqR, which is independent ofN. This is the same as EqUJfor Xo <^ 1. Thus, in this regime, 
ISLA predicts that mutational competition with the wild-type is the dominant factor limiting mu- 
tator fixation, and we recover the explicitly frequency dependent heuristic picture. In the opposite 
extreme of strong-effect mutators, regardless of the sign of a^iNs + 1) — 1, we recover our /i_ = 
result (Eqs Jl 11121) where deleterious mutations are the dominant factor limiting mutator fixation. 

These are pleasing mathematical results that seem to reconcile opposing heuristic viewpoints. 
However, they do not always match simulations in the weak-effect mutator regime. Figl6] (right) 
shows numerically generated solutions to Eql6] (Eq[T3] gives the large fi limit of these curves) as 
compared to the outcome of simulations. The disagreement is obvious: ISLA drastically overesti- 




(13) 



24 



mates the effect of the mutations in wild-type backgrounds. Figj6] shows that beneficial mutations 
in wild-type backgrounds eventually decrease for large enough R, though the decrease here 
is smaller than what ISLA predicts. The small effect of these mutations persisted even when we 
used parameters such that the wild-type subpopulation generated mutations at a rate A^(l — Xo)/i- 
that was equal to or even greater than the corresponding rate NxojJi+ in the mutator subpopulation. 
Although we do not fully understand this discrepancy, we can point to its source: There is a subtle 
error involving the final term of both Eqs J4|5| which states that during a single time-step, the mutator 
has a probability (1 — x)/i_aeS of becoming instantly lost. This is incorrect. The correct statement 
is that (1 — x);U„aeS is the probability that during one time-step, the wild-type generates a beneficial 
mutation that will eventually escape loss to random drift. Such mutations sweep through the popula- 
tion during a mean time interval t sweep ~ i^^^ generations which is typically much longer than the 
time to extinction of a mutator due to random drift tdrift ~ In(A^) (|Crow. J.F. AND KiMURA. M 



1970|) . However, for sufficiently large s, t sweep is small, Al becomes a better approximation, and 



ISLA more closely matches simulations. An example of this agreement is presented in Supple- 
mentary Information, where s = 1/3, = 1000, = .4, i? = 10. Thus, ISLA provides accurate 
results except in the weak-effect mutator regime with sufficiently small s. Unfortunately, we do not 
have a quantitative sense as to how large s must be in order to achieve accuracy. We plan to address 
this issue in future work. 

In Supplementary Information, we more closely examine the role of yU_ by presenting and in- 
terpreting the distribution of fixation and loss times for mutators when = and = 100. 



COMPARISON OF ISLA TO SIMULATION 

We now return to the case /i_ = 0, where the results of ISLA agree with simulations when 
R = n^/n^ is sufficiently large. Figs J4l5l illustrate the agreement between ISLA Eql6]and simu- 
lations, whenever is not too large. However, for larger /i+/s, we see the emergence of two 
qualitatively distinct discrepancies between ISLA and simulations. For < 1, a relatively small 
difference accumulates, whereas when ^+/s reaches values of 0(1), a drastic difference emerges. 
In this section, we analyze the sources of these discrepancies. 
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Figure 7; The role of non-lethal deleterious mutations. We "turned off" deleterious mutations, both in 
simulations and in ISLA, by setting the deleterious mutation rate to zero and leaving the beneficial 
mutation rate unchanged (left). The difference between these results and the corresponding ones 
with deleterious mutations is plotted on the vertical axis on the left. For < 1, deleterious 

mutations have the same effects in ISLA Eql6] as in simulations (left). ISLA essentially treats 
deleterious mutations as lethal (A3), instead of merely having a selective disadvantage —s. We 
tested this approximation directly in simulations by varying the parameters a and 6 while holding 
the product a{l — 5) = cxe constant (right). Parameters are s = 1/120, = 5000, = and 
a = A,5 = (left only). 
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The broad reason that ISLA and simulation do not agree for all yU+ is simply that A1-A3 and 
the resulting transition probabilities are only an approximation of the complex stochastic process 
executed by the simulations. Indeed, strictly speaking, the simulation does not even undergo a 
Markov process with respect to the variables x, t: one must also consider the fitness distributions 
of the subpopulations in order to write down the exact transition probabilities. When viewed this 
way, it is perhaps surprising that A1-A3 work as well as they do. We now specifically point out the 
errors introduced as a result of A1-A3, all of which are associated with mutational processes. 

A3 is accurate when n+fs < 1: We first analyze the way that ISLA treats deleterious mutations, 
which includes both A3 (which treats all deleterious mutations as lethal) and Al (which does not 
allow deleterious mutations to arise in the course of fixation of an "evolved" clone). FigjT] (right) 
compares simulation results from two sets of parameters with identical beneficial mutation rates 
(ae/i+) but different allocations of lethal and deleterious mutations via a difference in the parameter 
6. The results are essentially identical as long as /i+/s < 1. This shows that as far as mutator 
fixation is concerned, mutations of effect —s can be considered lethal, i.e. A3 is accurate in this 
regime. 

Al is accurate when ^+/s < 1: Furthermore, we can test all the effects of deleterious mutations 
by removing them from both the simulations and ISLA: the deleterious mutation rate is set to zero 
whereas the advantageous mutation rate is left unchanged. The results of this case are presented 
in FigjT] (left). Predictably, increases monotonically with /i+ in this case (data not shown). To 
compare the effect of deleterious mutations in simulations against those same effects according to 
ISLA, Eqil we plot the difference AS*^ = S^^no-deieterious - S ^^deleterious between results with 
deleterious mutations "off" and those with deleterious mutations "on" in the two cases. We see 
in FigjT] (left) that AS^j from ISLA matches that from simulation until /i+/s 1. Also note that 
AS^ ^ for ^ .1, illustrating the negligible effect of deleterious mutations in this regime. 
Thus, both Al and A3 are accurate when yU+/s < L 

A2 fails when jj+fs < 1: Since Al and A3 remain valid in this regime, the mild discrepancy 
between simulations and ISLA must originate in A2, which handles beneficial mutations. Specifi- 
cally, the fraction (1 — s) of advantageous mutants that are lost to random drift are treated as neutral 
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mutators which can later give rise to beneficial mutants that may sweep through the population. In 
some sense, this overstates the potentia l of these mutants because, in fact, they are typically lost to 



random drift within a few generations (|Crow. J.F. and Kimura. M 



19701). There is no simple 



remedy for this deficiency in A2, but an alternative, which we denote A2*, is to immediately kill 
these advantageous mutants, thereby treating them equivalently to deleterious and lethal mutants. 
Whereas A2 overestimates in this regime, A2* underestimates it. Thus, the simulation data is 
bounded by the predictions of A2 and A2* when ^ 1. See Supplementary Information for a 
graphical comparison and further discussion of A2*. 

Al fails when ji+j s ~ 1.- We now turn to the large discrepancy between ISLA and simula- 
tions when is 0(1), as seen in FigH Roughly speaking, this occurs when the time-scales 
of (deleterious) mutation and selection become comparable. In this regime, members of an ex- 
panding "evolved" clone are "lost" due to deleterious mutations faster than they are "added" due 
to selection. Consequently, the fixation probability of an advantageous mutant in a homogeneous 
genetic background 7r(s) < s and Al fails. Semi-quantitatively, we expect this effect to set in when 
(1 — ae)iJi+/ s ~ 1. The dependence can be seen by comparing Figs[ 



COMPARISON TO EXPERIMENT 



As mentioned previously, the spontaneous emergence of mutator alleles has been documente d in 



laboratory evolution experiments with E.coli (|SnieG0WSKI et ali 



1997 



Shaver et al. . 



2002h . In 



this experiment, mutator alleles with R ^ 100 became fixed in 3 out of 12 independently evolving 
E.coli populations within 10,000 generations. The total number of mutators generated among 12 



lines during 10,000 generations is approximately Ne x U x (10^ x 



rate into the mutator state and is the effective populatio n size (IWahl and Gerrish 



Wahl et al. 



(BOEet al. 



2), where U is the mutation 



2001 



1997ah and 5 x 10' 



20021) . U has been measured between 5 x 10 (^Taddei et al. . 
20001) . and we find N^, = 6.3 x 10'' (see Supplementary Information). Since three of 



these mutators achieved fixation, the experimental fixation probability Pfix,expt is approximately 
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given by 3/(A^'e x U x 10^ x 12) and bounded by 



7.9 X 10-« < Pfi,^,,pt < 7.9 X 10-^ (14) 

This value is 5-50 times that of a neutral allele (1 /N^,). 

In order to compare this value to the predictions of ISLA, we need experimental values for the 
parameters ae, and s. It turns out that the equivalent set of parameters s, the beneficial mutation 
rate fj,ben,+ = ctefJ'+, and the deleterious mutation rate fidei,+ = (1 — «e)/^+ are more readily 
available in the literature. A survey of these parameter values, as well as a more careful discussion 
of their meaning, can be found in Supplementary Information. Presently, we use the beneficia l 



mutation r ate Uhen = 2.8 x 10 and selection coeffi cient s = .1 obtained by 



Lenski et al. 



Following 



mil) 



Keightley and Eyre-WalkerI (|1999|) . we take fidei = 1-6 x 10" ^ These mutation 



rates are based on the measured wild-type values and assume R = 100. Since Nefidei,+ ^ 1, 
N^lj,hen,+s ^ 1, and N^OeS <^ R, these populations are in the drift-less, strong-effect mutator 
regime. Therefore, the appropriate formula is either Eq[T2lor Eq[T3l which give the same results. 
Plugging our parameter values into ISLA, we obtain 

PfixMa = 1-8 X 10-« (15) 

in reasonable agreement with the rough experimental value (Eq [T4)) . Other choices for parame- 
ter values, particularly fiben,+, would result in less impressive agreement with experiment. See 
Supplementary Information for further discussion. 

It is also interesting to note that, according to these experimental parameters, Na^s ~ 1.1, 
indicating that these E. coli populations only very marginally favored mutators. This could ex- 
plain why no mutators fixed during the next 25, 000 generations: NueS had decreased below the 
threshold value of one as fewer, and less potent, beneficial mutations became available. 

Due to the relatively large population size = 6.3 x 10^ and the anticipated small fixation 
probability, we cannot obtain an accurate measurement of Pfix using our simulation method. How- 
ever, for these experimental parameters, /i+(l — ae)/s = fidei,+/s is 0(1) and therefore we expect 
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the data to lie in the decreasing portion of curves such as FigjSl Thus, our ISLA estimate of Pji^ 
is probably much larger than what simulations would yield. We briefly return to this issue in the 
Discussion. 



DISCUSSION 

Relation to Previous Theoretical Work: As mentioned in the introduction, there are many 
existing theoretical models of mutator evolution. In this section we briefly review the existing body 
of knowledge and place our present work in this larger context. Studies are discussed roughly in 
order of increasing similarity to our present w ork. 



Models with explicit environmental change: 



Leigh JrI (119701) endeavored to calculate the muta- 



tion rate that maximizes the growth rate of its corresponding modifier locus. An infinite population 
with this wild-type ("resident") mutation rate is evolutionarily stable in the sense that it cannot be 
invaded and swept by any modifier of mutation rate. Su ch an evolutionarily stable strategy (ESS) 



is referred to as the ESS mutation rate. 



Leigh JrI (|1970|) developed a simple two locus, two allele 



model of mutator dynamics in an environment that regularly alternates between two states. One lo- 
cus is under selection, and its two alleles are alternately favored in the two different environments. 
The second locus is not under direct selection and merely modifies the mutation rate at the selec- 
tive locus. The dynamics of the mutator allele are deterministically governed by two effects. First, 
immediately after the environment changes, the mutator increases its frequency because the small 
population of mutants, which is favored in the new environment, is over-represented in the mutator 
background. This favors the higher mutation rate. Secondly, after the mutant sweeps through the 
population, the frequency of the mutator decreases due to association with the deleterious mutants 
that it generates at its new fitness peak. This favors lower mutation rate. The cycle repeats itself 



many times, and 



Leigh JrI (|l970i) finds that the long term ESS mutation rate is equal to the rate of 



environmental change. Over the years, this basic model was improv ed by incorporatin g the effects 



of timing of en vironmental changes 



varyin g selective coefficients ( 



ate ge notypes (|Travis and Travis 
2006I) . 



ISHH et al. 



2002h . and multiple mutable sites (|Palmer and LlPSlTCH . 



1989). intermedi- 
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Figure 8: The scaling behavior of EqU] and ISLA are qualitatively different. If the initial number 
of mutators Nxo is kept constant while N is increased, then ISLA predicts that Pjj^. remains in- 
variant, whereas the frequency dependent Eq{T] predicts a large change. Simulations are in better 
accord with ISLA than Eq{Tl These scaling predictions could be experimentally tested by observ- 
ing whether the "threshold" number of initial mutators changes with N . Here, we have defined the 
threshold as the number of mutators for which Pjix = 1/2, and depicted these values with vertical 
dotted lines. Parameters are a = .4, 5 = 0, /i+ = s = 1/120. 
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While these models doubtless provide valuable insight into certain biological scenarios, they are 
rather orthogonal to our work. Three differences seem especially important. First, most obviously, 
mutator success requires repeated environmental changes in these models. In contrast, our model 
shows that environmental change is only necessary for mutator fixation insofar that it provides a 
rationale for having a population displaced from its fitness peak. Secondly, they endeavor to find 
the global ESS mutation rate whereas we focus on quantifying, via fixation probability, the proba- 
bilistic result of a single competition experiment. While full knowledge of Pfi^^N, s, a, S) 
implies the value of the ESS, the converse is not true. Thirdly, their mechanism of mutator success 
is very different from ours. Whereas they rely upon the alternating selective effects of existing 
mutants to boost mutator frequency, our model analyzes the dynamic, stochastic interplay between 
random drift, deleterious mutations, and advantageous mutations in a constant environment. We 
propose that, on the whole, our model contains fewer special assumptions than models with ex- 
plicit environmental change. Regardless of whether fluctuating or constant environments are more 
biologically informative, our results consti tute an important null model of mutator fixation. 



Constant environment models: Work by 



Tanaka et al. 



(|2003|) also involves a changing environ- 



ment. However, unlike the models described in the previous section, theirs contains no alternating 
selective effects: when the environment changes, the mutations acquired during the previous envi- 
ronmental cycle simply become neutral. Thus, as in our work, all ben eficial mutants are generated 



de novo. In further similarity with our work. 



Tanaka et al. 



(120031) pursue, via quasi-stochastic 



simulations and analytic approximations, an understanding of the long term mutator behavior by 
concentrating on a single environmental cycle, i.e. by examining populations in a constant environ- 
ment. These authors were interested primarily in the case when A^Xo/i+ -C A^(l — Xo)/i-, where the 
fixation of mutators is in some sense unlikely. With this in mind, instead of Pfix, they measure and 
calculate the (much larger) probability Pgain that the initially rare mutator increases its frequency by 
the end of a "time cycle." These cycles are defined to end when an expanding clone in a wild-type 
background reaches a size of 0{N), at which point the simulation is halted. Their most interesting 
result is that Pgain is substantial even when NxoiJ.+ <^ A^(l — Xo)iJi^. In other words, mutators 
can still "break even" if the wild-type background generates the first beneficial mutation, which is 



32 



important if the environment changes. Nonetheless, without environmental change in their model, 
mutators will always be doomed unless they are the first to generate a beneficial mutation. Further- 
more, they model birth and death processes deterministically, in a way that precludes extinction. 
For these reasons, our Pfix and their Pgain are truly distinct quantities, and no direct comparison 
can be made with our work. 

We next discuss a simple calculation by 



Lenski 



( 2004 ) based on indirect mutation- selection 



equilibrium of the mutator subpopulation. If the dominant processes occurring in the population 
are mutation into the mutator state and creation of deleterious mutations by mutators, then the 
frequency of mutators approaches an equilibrium value. This frequency is easily calculated if, as 
in A3 of ISLA, deleterious mutations are treated as immediately lethal: 



U 



U 



X, 



eq 



The time taken for the population to reach this equilibri um state, as well as a much more careful 



calculation of Xeq, was investigated by 



Johnson! (|1999bh . but presently we assume that this simple 



estimate is sufficient. In equilibrium, beneficial mutations therefore arise at a rate NxeqjJi+cie from 
the mutators, and rate A^(l — Xeg)jj^ae from the wild-type. If all beneficial mutants of equal effect 
have the same probability of achieving fixation, regardless of whether they originate in a mutator 
or wild- type background, then the fraction of substitutions linked to a mutator is approximately 



U 



u 



(16) 



Plugging in reasonable values. 



LenskiI (20041) finds that 1% of substitutions should be linked to 



Sniegowski et al. 



(119971) 



mutators. Furthermore, given that each line of E. coli in experiments by 
generated 10-20 substitutions, this calculation is impressively consistent with the observation that 
3/12 lines became mutators. 

In order to relate this approach to our own, we must reintroduce dynamics into the picture. We 
can interpret the quantity as the conditional probability that a mutator achieves fixation, given 
that a selective sweep occurs during its lifetime. Our quantity Pfi.j. is this conditional probability 
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multiplied by the probability that a selective sweep occurs during the lifetime of a mutator. As- 
suming that selective sweeps and death each occur as Poisson processes with rates NfxaeS and 
— — «e), respectively, it is straightforward to show that the probability that at least one 

selective sweep occurs before death is given by 



(1 - ae)(/i+ - /i-) V (1 - «e)(/i+ - A^-) 

Multiplying this expre ssion by the con ditional probability ~ xR, we obtain Eq{T3l Thus, the 



r 

approach suggested by iLenski 



(|2004l) is the equilibrium version of ISLA, in the limit where mu- 
tational processes occur frequently enough to overwhelm random genetic drift. Thus, remarkably, 
even though this approach frames the problem of mutator fixation in terms of competition with 
beneficial mutations in wild-type backgrounds, R cancels out of the solution in the strong-effect 
mutator regime: R ^ Naes/{1 — Oe). 

It is also worthwhile to examine the conditions under which we expect the equilibrium as- 
sumption to hold. Let us imagine that an evolution experiment is conducted for T generations, 
during which H substitutions occur. ISLA predicts that the expected number of mutator fixations is 
NPfiJJT, whereas according to EqdH the equilibrium approach yields a value equal to Hj. Set- 



ting these two values equal to one another, and plugging in (from Eq[TT]) Pfi^ {xo = 1/N) = Sr^^, 
we obtain 

H = NsfiT ^ NfiaesT 

1 - tte 

This expression merely states that the (mostly wild-type) population is in the "successive mutations 
regime", i.e. only a single beneficial mutation spreads at a time. Alternatively, one could imagine 
turning this argument around and asking what Pfi^ must equal given that the equilibrium approach 
is valid and that the population accumulates substitutions "one by one". In that case, one would, 
remarkably, arrive at Pfix{xo = ^/N) = aeS, which (for small ae and NS^ ^ 1) is what we 
obtained earlier (Eq fTT]) by r nore sophisticated method s. 



Turning to another study. 



Tenaillon et al. 



(|1999|) investigated, via stochastic simulations and 
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very brief analytic arguments, multi-locus mutator evolutio n in a constant enviro nment. These 



Taddei et al. 



(1 19970) and are partly 



extensive simulations are a generalization of earlier work by 
amenable to comparison with our work. Some noteworthy differences with our simulations are 
that they scan a larger range of N, they have a more realistic implementation of mutation, and, 
most importantly, they allow flux into and out of the mutator state. Thus, mutators are never 
absolutely fixed during their trials, which necessitates a different termination condition than ours: 
They declare a trial "over" when the population reaches its maximum fitness, whereas we declare it 
"over" when the mutator is completely and permanently fixed or lost. Upon termination of the trial, 
they consider the mutator "fixed" if its frequency is > 95%. They measure the fraction of trials 
that terminate with mutator frequency > 95% and denote this quantity the "frequency of mutator 
fixation," which differs from our Pfi^ because of reasons discussed below. 

One important consequence of their method is that the total number of mutators generated dur- 
ing a trial varies with the choice of parameters. This is because each replication event presents a 
chance for the creation of a new mutator, and the number of replication events that occur before 
termination clearly depends on N, s, and the number of mutational steps required to reach 

the peak. Thus, a change in the value of any of these parameters may alter the "frequency of mu- 
tator fixation" simply because it changes the number of mutators that are typically created during 
the trial. Our Pfix, on the other hand, remains invariant under such changes and allows us to filter 
out this background effect. Their system is doubtless a more literally accurate representation of 
biological reality, which has its virtues but also major costs, which we discuss below in the context 
of two important examples. 

First, they measure that the "frequency of mutator fixation" increases with A^. This is an in- 
teresting and potentially practical result, but their method makes it very difficult to determine the 
extent to which the increase is simply due the background effect that more mutators were created in 
the larger populations. ISLA, on the other hand, unambiguously states that when NS^ ^ 1, Pfix 
for a single mutator becomes independent of N. Therefore, ISLA predicts that the d ependence 



of mutator fixation frequency on population size observed by 
driven by the simple background effect. 



Tenaillon et al. 



(| 19991) is entirely 
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A second example has even more dramatic conceptual consequences. These authors ask whether 
Pfix is determined by the number of potentially advantageous mutations (steps away from the peak) 
or merely by the rate that such mutations are generated. In order to investigate this question, they 
devised two sets of simulations. In one set, there were 12 available advantageous mutations, ac- 
cessible at a rate of 10"*^ each. In the other set, there was a single mutation of the same effect, 
accessible at a rate of 12 x 10~®. The explicit difference between these sets of simulations is the 
number of steps to the fitness peak, but an additional, implicit difference is that the set with 12 
beneficial mutations runs for more generations. Therefore, more mutators are created in that set of 
simulations. Now, ISLA predicts that Pfix depends only on the advantageous mutation rate, and 
that therefore the two simulations should result in the same Pfix- In seeming contrast, they found 
the "frequency of mutator fixation" to equal approximately .5 for the first situation and approxi- 
mately zero for the second. This observation led them to conclude that mutators succeed because 
of their advantage in rapidly creating genomes which carry multiple beneficial mutations, which 
is fundamentally different from our conceptual picture. We propose that this simulation finding 
might be explained by the simple background effect that far more mutators are created en route to 
acquiring 12 beneficial mutations than to acquiring a single beneficial mutation. ISLA completely 
neglects multiple beneficial mutations, and its success, both nea r the peak (Fig]4l) and far from it 



Tenaillon et al. 



(1999) in fact 



(Figj5]), suggests that the multiple mutations effect proposed by 
plays a very minor role in mutator fixation. However, it should be noted that we did not investigate 
cases where the mutator is favored and only a single beneficial mutant is available. It could be the 
case that multiple beneficial mutations in the same genome are implicitly important in that they 
are what allows the mutator to overcome competition with wild-type beneficial mutations. This 



lould be explored in futur e work 



hypothesis s. 

Whereas 
by I Andre and Godelle 



Tenaillon 



et al. 



many similarities to ours. 



1999|) focused almost exclusively on stochastic simulations, work 



(120061) relies almost exclusive ly on analytic methods. In work that bears 



Andre and GodelleI (|2006f) studied, mostly via an analytic approach. 



the long term trajectory of mutation rate evolution. A key insight of theirs is that, in a finite asexual 
population, the frequency of a mutator undergoes strong fluctuations, with values covering the 
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entire range from zero initially to one upon a selective sweep by a linked locus. Thus, they point out 



that studies which assume that mutators are rare during all generations, either because of infinite 
I 1 I 1 I ■ ■ ■ 



population size (iLeigh Jri . 1 19701) or sex ual recombination (iJOHNSON 



different than finite asexual populations. 



Andre and GodelleI (|2006|) remedy this problem by 



1999a|) . are qualitatively 



calculating the fixation probability of an initially rare mutator. We now briefly summarize their 
method of solution and show that, with minor modification, it corresponds to the Nfi oo limit of 
our results. In what follows, we take some liberty in changing their notation and using continuous 
time. 

Their initial condition is identical to ours: a clonal population is seeded with a small number of 
otherwise identical mutators. They then temporarily ignore beneficial m utations and analyz e how 



the expected number of mutators changes with time. In agreement with 



JohnsonI (|1999bl) . they 



find that after a waiting time 1/s, the mutator subpopulation declines exponentially, i.e. E[x(t)] 



~{^l+~^l.){l~ae){t~l/s) 



. They then construct their key equations (their Eq. 19) 



d 



;i - Pf,,{t) - Piossit)) ■ NfiaeS ■ ^ ■ E[x{t)] 

;i - Pfi,{t) - Ploss{t)) ■ NfxaeS ■ ^ ■ (1 - E[x{t)]) 



We have written these equations in a somewhat peculiar way, and replaced their symbol K with 
NftaeS in order to facilitate translating between our notation and theirs. These equations are very 
similar to ISLA in that they represent the instantaneous fixation of beneficial mutations which orig- 
inate from a time dependent mutator subpopulation. However, there are two disturbing features 
about these equations. First, they assume that the only cause of mutator extinction is beneficial 
mutations in the wild-type background. In fact, mutators also become extinct due to (i) their muta- 
tional load and (ii) random drift. In their equations, E[x{t)] declines exponentially, but erroneously, 
this decline does not contribute to Pioss- Both (i) and (ii) cause an overestimate of Pfix- The second 
distur bing feature of these equations is the appearance of expectation values on the RHS. With this 



move. 



Andre and GodelleI (|2006r) replaced the random variable x{t) with its mean value, which 



is a very substantive approximation. The distribution of x{t) is in fact diffusing, i.e. random drift is 
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in fact occurring. Nevertheless, we expect that their representation of x{t) as a deterministic quan- 
tity to be approximately valid when the timescale of this diffusion is slower than the timescales due 
to mutation and selection. Unlike our approach, theirs cannot quantify when it is safe to neglect 
random drift. Looking back to Eql6l we see that the diffusive process, i.e. random drift, can be 
neglected when Nfi±{l — a^) ^ 1 and N'^^^a^s ^ 1. It just so happens that these criteria will 
often be met in microbial populations. 



We now explicitly demonstrate some important parallels between our work and that of 
(|2006|) in the large Nfj, limit. Since, in our model, deleterious mutations are as strong as advanta- 
geous ones, the best comparison is made with their "ruby in the rubbish" hypothesis. The relevant 
solution is their Eq.A5 



Andre and Godelle 



Pfi 



fix 



NjjiaeS 



(17) 



Simplifying the denominator by taking exp [—{^+ — — Og)] ~ 1 — (/Lt+ — — ote) and 

neglecting the term —NjlaeSH+, we recover our large Njx result from ISLA (EqfT3l). The neglected 
term inflates the value of Pf ix, and is a result of these authors not treating extinction of the mutator 
due to its mutational load. This has important consequences for the next topic. 

Long term mutation rate evolution: Although our work primarily addresses the plain issue 
of calculating Pji^, we briefly contemplate implications for the more grand question of long term 
mutation rate evolution. 

^J'conv is proportional to the rate of sweeps: Thus far we have considered selective sweeps to 
be initiated by de novo beneficial mutations. Let us now briefly apply our results to the case 
where sweeps are instead triggered by an environment that changes at ra te K. This merely re- 



quires transcribing NjiaeS ^ K. Following 



Andre and GodelleI (l2006h we expand the fixation 



probability (EqilTl) in powers of /i+ — fi^ and denote the first order coefficient in this series by 
Sel{jj-.). The roots of Sel{n-) give the "convergence stable resid ent mutation rate ." Using EqfT3l 
we find Urnr,,, = K/(l — ftp) ^ K, which is the classical result (|Leigh JrI . I1970|) . Using EqfTTl 
Andre and GodelleI (|2006|) find a qualitatively different result: Umn,, = -^^^ n j^y w hich di- 



verges as K 1. The reason for this discrepancy is that 



Andre and Godelle 



( 2006 ) did not 
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allow for extinction due to the mutational load. ISLA naturally accounts for these extinction events 
and leads to the classical result. However, ISLA approximates deleterious mutations as being lethal, 
whereas these authors also treated the more realistic non-lethal case. It may be possible to demon- 
strate, via further analysis, the claim that non-lethal deleterious mutants cause Hconv to diverge for 
some parameter values. 



Equilibrium mutation rate: We find that Selj^iJ] = ~~ ^^a^J ^ ' whereas 



Andre and Godelle 



(|2006h find Sel{fi-) = J_M^±1)_J. _|_ i _ j^^. Our expression indicates that there are no equilib- 



rium mutation rates: for all weak mutators are favored when ae{Ns + 1) ^ NaeS > 1 and 
disfavored in the opposite case. This threshold is clearly in agreement with our EqjTl Thus, as far 
as ISLA is concerned, populations with NaeS < 1 should continually evolve toward the minimum 
attainable mutation rate. On the other hand, populations with NagS > 1 should evolve an ever 
higher mutation rate. Our expression for Sel{fj._) is clearly inaccurate for very small /i_ (because 
random drift dominates in that regime) and also for very large /i_ (since our simulations show that 
there is a maximum mutation rate that can achieve fixation). 

Limitations of Present Work: Real biological populations possess many features that this 
article either neglects or severely constrains. We now briefly discuss the most striking limitations. 

Initial Conditions: Both ISLA and our simulations suppose that "initially" all members of the 
population have the same fitness. If this assumption is false and mutators arise randomly in a 
population with pre-existing fitness variation, this might act to decrease mutator success: unless the 
mutator happens to emerge from the fittest subclass of the population, the advantageous mutations 
it generates will already be present in more abundant subclasses which could out-compete the rare 
mutator. This point is especially relevant since, in comparing ISLA to experiment, we essentially 
assumed that each mutator that arose during the course of the experiments did so in a population 
consisting of a single fitness value. 

Strict Asexuality: Our simulations and ISLA do not allow any mechanisms of horizontal gene 
transfer or recombination. These events would decouple mutator alleles from the advantageous 
mutations that they generated, and thereby result in significantly decreased mutator success. This 
effect is especially important since some genes associated with a mutator phenotype also exhibit 
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hyper-recombination (|Denamur and MATld . 



2006 ). 



Simple Fitness Landscape: Our simulations assume that mutations all fall into one of three 
classes: lethal, beneficial with effect +s, or deleterious with effect —s. As mentioned previously, 
and discussed in Supplementary Information, it may be true that, in large pop ulations, beneficial 



mutations of a fixed size s are the ones that typ i cally r each appreciable frequency (IGerrish and Lenskii . 



1998 



Desai et al. 



2007 



Hegreness et al. 



20061). However, this simplification is certainly not 



possible when considering deleterious mutants, whose distribution is 



modal, with many mutations being nearly neutral and many being lethal (|Eyre-Walker and Keightley 



ikely complicated and bi- 



20071) . Fig j7] suggests that increasing the strength of deleterious mutations has effects only at large 



yU+/s, where it increases both the peak value of and the value at which the peak oc- 

curs. Along these lines, a simulation model that included a class of weakly deleterious mutations 
would likely continue this trend. This would delay the large discrepancy between the simulations 
and ISLA until even large r u^l s. This issue could h elp to explain the previously mentioned fact 
that in experiments of ISniegowski et al\ (|l997h seem very close to the maximum allowable 
value. Including mildly deleterious mutations would also prolong the lifetime of genomes which 
carry them. In this case, it might be necessary to inc orporate a time dela y before these deleterious 



Johnson! (1999b). 



mutations are "enforced," along the lines explored by 

Suggestions for further research: This article leaves many questions unanswered, but also 
points to interesting theoretical and experimental opportunities. 

Theoretical directions: A satisfactory analytic description of our stochastic simulations remains 
incomplete. Two key issues remain unresolved. First, we do not understand the mechanism by 
which mutators continue to succeed when faced wi th intense mutational competiti on from the wild- 



type background (Fig©. Our work and that of 



Andre and GodelleI (|2006h both imply that 



mutations in wild-type backgrounds should become important when NaeS ~ jj,-, but this is 
not borne out in the simulations unless s is "sufficiently large." Secondly, it is clear that ISLA fails 
to match simulations when the mutation rate is very large (1 — ae)/i+ ^ s. Quantifying the success 
of mutators in this regime is especially relevant to studies of long term mutation rate evolution. 
Another issue that we did not address is the full dynamics of mutator fixation. Our analytic 
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results are mostly derived from Eql6l which is relevant to the eventual fate of mutators. An ap- 
proximate solution to the time dependent forward diffusion (Eqj?]), with /i_ = 0, is given in Sup- 
plementary Information. This solution provides some dynamical information, but, like the entire 
ISLA approach, it assumes that selective sweeps occur instantaneously. In this sense, Eq|4]predicts 
incorrect dynamics. Furthermore, we showed that mutator success is compactly represented by an 
effective selection coefficient S^j. For simple advantageous mutants, S contains information not 
only about Pfi^ but also about the average dynamics: ~ e'^* when rare. Perhaps that is the 

case with mutators as well. 

Experimental ideas: Our work shows that, in most regimes, Pfi^ is not explicitly frequency 
dependent. Rather, Pjj^. depends on the initial number of mutants Nxo- This scaling behavior 
could be tested experimentally. Suppose that competition experiments in a chemostat carrying a 
population of size A^^i showed that, when the initial frequency of mutators exceeded a threshold 
value of xi, mutator achieved fixation with a high probability. One could decrease the population 
size to N2 and again inoculate with mutators at a frequency of xi. Our results predict that mutators 
would not achieve fixation in this case because N2X1 is le ss than the threshold n umber A'^iXi. In 
fact, very similar experiments were recently performed by ILe Chat et al\ (|2006h . which support 
the notion that Pfix scales with Nxo and not with Xo alone. However, these competition experiments 
were done under a lethal selective pressure, which selected for pre-exiting resistant mutants. Here 
we propose competitions between initially isogenic (aside from the mutator allele) mutator and 
wild-type strains adapting to a new environment. In addition to this scaling behavior, ISLA predicts 
a testable value for this threshold that differs significantly from the frequency dependent picture 
represented by EqU] These ideas are presented in FigjSl 

It would also be interesting to experimentally investigate the decline in mutator success seen 
for very large mutation rates when (1 — ap )u+ ~ s. As ment i oned p reviouslv, during the first 



few thousand generations of experiments by 



Sniegowski et al. 



(Il997h . NaeS ^ 1.1. The reason 



why no mutators achieved fixation after the first 10, 000 generations could be that this parameter 



decreased below the thresho 
was previously discussed by 



d value of one during the co urse of its evolution. A similar effect 



Kessler and LevineI (|1998|) . An alternative explanation is that ;U+ 
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wa s near the theoretical maximum (1 — ae)/^+ ~ 1 suggested from our simulations. As noted 



by 



Gerrish et al. 



(120071) . once could test these competing explanations by founding several new 
lineages with a clone from of one of the mutator populations, and growing these mutator lineages 
in a novel environment. The new environment should be one in which NapS > 1. If no "double 



mutators" arose, then the hypothesis of a maximum allowable mutation rate would be supported. 
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